Chihik's Blog

An Ac a day, keeps the doctor away!

OJ6388 gcd求和

发布于 2020-10-14 | 分类于 莫比乌斯反演 | 13分钟 | 1982字数 | 浏览量：

data 1~10

令s=min(n,m)\text{令s=min(n,m)}

∑i=1n∑j=1mgcd(i,j)\sum_{i=1}^n\sum_{j=1}^mgcd(i,j)

∑k=1sk∑i=1n∑j=1m[gcd(i,j)=k]\sum_{k=1}^sk\sum_{i=1}^n\sum_{j=1}^m[gcd(i,j)=k]

∑k=1sk∑i=1⌊nk⌋∑j=1⌊mk⌋[gcd(i,j)=1]\sum_{k=1}^sk\sum_{i=1}^{\lfloor \frac{n}{k} \rfloor}\sum_{j=1}^{\lfloor \frac{m}{k} \rfloor}[gcd(i,j)=1]

∑k=1sk∑i=1⌊nk⌋∑j=1⌊mk⌋∑d∣gcd(i,j)μ(d)\sum_{k=1}^sk\sum_{i=1}^{\lfloor \frac{n}{k} \rfloor}\sum_{j=1}^{\lfloor \frac{m}{k} \rfloor}\sum_{d|gcd(i,j)}\mu(d)

∑k=1sk∑d=1⌊sk⌋μ(d)⌊nkd⌋⌊mkd⌋\sum_{k=1}^sk\sum_{d=1}^{\lfloor \frac{s}{k} \rfloor } \mu(d) \lfloor \frac{n}{kd} \rfloor\lfloor \frac{m}{kd} \rfloor

这样可以用两层分块 Θ(n)\Theta(n) 求 ， 复杂度Θ(n+qn)\Theta(n+qn)

#include <cstdio>
#include <iostream>
using namespace std;

const int MAXN = 1000000;
int t , n , m , k , x , prime[ MAXN + 5 ] , mu[ MAXN + 5 ] , f[ MAXN + 5 ];
bool vis[ MAXN + 5 ];

void sieve( ) {
	mu[ 1 ] = 1;
	for( int i = 2 ; i <= MAXN ; i ++ ) {
		if( !vis[ i ] ) {
			prime[ ++ k ] = i;
			mu[ i ] = -1;
		}
		for( int j = 1 ; j <= k && 1ll * i * prime[ j ] <= MAXN ; j ++ ) {
			vis[ i * prime[ j ] ] = 1;
			if( i % prime[ j ] == 0 ) break;
			mu[ i * prime[ j ] ] = -mu[ i ];
		}
	}
	for( int i = 1 ; i <= MAXN ; i ++ )
		f[ i ] = f[ i - 1 ] + mu[ i ];
}

long long Get( int n , int m ) {
	long long Ans = 0;
	for( int l = 1 , r ; l <= min( n , m ) ; l = r + 1 ) {
		r = min( n / ( n / l ) , m / ( m / l ) );
		Ans += 1ll * ( f[ r ] - f[ l - 1 ] ) * ( n / l ) * ( m / l );
	}
	return Ans;
}
long long solve( int n , int m ) {
	long long Ans = 0;
	for( int l = 1 , r ; l <= min( n , m ) ; l = r + 1 ) {
		r = min( n / ( n / l ) , m / ( m / l ) );
		Ans += 1ll * ( r - l + 1 ) * ( l + r ) / 2 * Get( n / l , m / l );
	}
	return Ans;
}

int main( ) {
	sieve( );
	scanf("%d %d",&x,&t);
	while( t -- ) {
		scanf("%d %d",&n,&m);
		printf("%lld\n", solve( n , m ) );
	}
	return 0;
}

data 11 ~ 15,21 ~ 25

这个是 part1\text{part1} 的优化，

∑k=1sk∑d=1⌊sk⌋μ(d)⌊nkd⌋⌊mkd⌋\sum_{k=1}^sk\sum_{d=1}^{\lfloor \frac{s}{k} \rfloor } \mu(d) \lfloor \frac{n}{kd} \rfloor\lfloor \frac{m}{kd} \rfloor

∑k=1sk∑kd=1sμ(d)⌊nkd⌋⌊mkd⌋\sum_{k=1}^sk\sum_{kd=1}^s \mu(d)\lfloor \frac{n}{kd} \rfloor\lfloor \frac{m}{kd} \rfloor

令T=kd\text{令T=kd}，

∑k=1sk∑k∣Tsμ(Tk)⌊nT⌋⌊mT⌋\sum_{k=1}^sk\sum_{k|T}^s \mu(\frac{T}{k}) \lfloor \frac{n}{T} \rfloor\lfloor \frac{m}{T} \rfloor

∑k=1s∑k∣Tsk×μ(Tk)⌊nT⌋⌊mT⌋\sum_{k=1}^s\sum_{k|T}^s k \times \mu(\frac{T}{k}) \lfloor \frac{n}{T} \rfloor\lfloor \frac{m}{T} \rfloor

∑T=1s∑k∣Tk×μ(Tk)⌊nT⌋⌊mT⌋\sum_{T=1}^s\sum_{k|T} k \times \mu(\frac{T}{k}) \lfloor \frac{n}{T} \rfloor\lfloor \frac{m}{T} \rfloor

∑T=1s(∑k∣Tk×μ(Tk))⌊nT⌋⌊mT⌋\sum_{T=1}^s(\sum_{k|T} k \times \mu(\frac{T}{k})) \lfloor \frac{n}{T} \rfloor\lfloor \frac{m}{T} \rfloor

括号里面的是卷积形式，id∗μ=φid*\mu=\varphi , 即

∑T=1sφ(T)⌊nT⌋⌊mT⌋\sum_{T=1}^s \varphi(T) \lfloor \frac{n}{T} \rfloor\lfloor \frac{m}{T} \rfloor

现在就可以一次数论分块解决了, 复杂度Θ(n+qn)\Theta(n+q\sqrt n)

#include <cstdio>
#include <iostream>
using namespace std;
#define Int register int

const int MAXN = 10000000;
int x , t , n , m , k , prime[ MAXN + 5 ] , phi[ MAXN + 5 ];
long long f[ MAXN + 5 ];
bool vis[ MAXN + 5 ];

template<typename _T>
void Read( _T &x ) {
	x = 0; int f = 1;
	char s = getchar( );
	for( ; s < '0' || s > '9' ; s = getchar( ) ) f = s == '-' ? -f : f;
	for( ; s >= '0' && s <= '9' ; s = getchar( ) ) x = ( x << 3 ) + ( x << 1 ) + s - '0';
	x *= f;
}
template<typename _T>
void Write( _T x ) {
	if( x < 0 ) putchar('-') , x = -x;
	if( x >= 10 ) Write( x / 10 );
	putchar( x % 10 + '0' );
} 

void sieve( ) {
	phi[ 1 ] = 1;
	for( Int i = 2 ; i <= MAXN ; i ++ ) {
		if( !vis[ i ] ) {
			prime[ ++ k ] = i;
			phi[ i ] = i - 1;
		}
		for( Int j = 1 ; j <= k && 1ll * i * prime[ j ] <= MAXN ; j ++ ) {
			vis[ i * prime[ j ] ] = 1;
			if( i % prime[ j ] == 0 ) {
				phi[ i * prime[ j ] ] = phi[ i ] * prime[ j ];
				break;
			}
			else
				phi[ i * prime[ j ] ] = phi[ i ] * ( prime[ j ] - 1 );
		}
	}
	for( Int i = 1 ; i <= MAXN ; i ++ )
		f[ i ] = f[ i - 1 ] + phi[ i ];
}

long long solve( int n , int m ) {
    int d = min( n , m ); long long Ans = 0;
	for( Int l = 1 , r ; l <= d ; l = r + 1 ) {
		r = min( n / ( n / l ) , m / ( m / l ) );
		Ans += ( f[ r ] - f[ l - 1 ] ) * ( n / l ) * ( m / l );
	}
	return Ans;
}

int main( ) {
	sieve( ); Read( x ) , Read( t );
    while( t -- ) {
        Read( n ) , Read( m );
        Write( solve( n , m ) ) , putchar( '\n' );
    }
	return 0;
}

data 16~20

可以参考一下P5739 LJS的GCD

∑i=1n∑j=1ngcd(i,j)\sum_{i=1}^n\sum_{j=1}^ngcd(i,j)

2∗∑i=1n∑j=1igcd(i,j)−∑i=1ngcd(i,i)2*\sum_{i=1}^n\sum_{j=1}^igcd(i,j)-\sum_{i=1}^ngcd(i,i)

2∗∑i=1n∑j=1igcd(i,j)−n∗(n+1)22*\sum_{i=1}^n\sum_{j=1}^igcd(i,j)-\frac{n*(n+1)}{2}

令 f(n)=∑i=1n∑j=1igcd(i,j)f(n)=\sum_{i=1}^n\sum_{j=1}^igcd(i,j)

f(n)=∑i=1n∑j=1igcd(i,j)f(n)=\sum_{i=1}^n\sum_{j=1}^igcd(i,j)

f(n)=∑i=1n∑d∣id∑j=1i[gcd(i,j)=d]f(n)=\sum_{i=1}^n \sum_{d|i}d \sum_{j=1}^i [gcd(i,j)=d]

f(n)=∑i=1n∑d∣id∑j=1id[gcd(id,j)=1]f(n)=\sum_{i=1}^n \sum_{d|i}d \sum_{j=1}^{\frac{i}{d}} [gcd(\frac{i}{d},j)=1]

f(n)=∑i=1n∑d∣id×φ(id)f(n)=\sum_{i=1}^n \sum_{d|i}d \times \varphi(\frac{i}{d})

熟悉卷积的可以发现，∑d∣id×φ(id)\sum_{d|i}d\times\varphi(\frac{i}{d}) 是一个标准的卷积式。

因为 id,φid , \varphi 为积性函数 ，所以卷积也为积性函数，考虑线筛

f(p)=2p−1f(p)=2p-1

f(pk)=∑i=0kφ(pi)×pk−if(p^k)=\sum_{i=0}^k \varphi(p^i) \times p^{k-i}

                     =pk+∑i=1k(pi−pi−1)×pk−i~~~~~~~~~~~~~~~~~~~~~=p^k+\sum_{i=1}^k (p^i-p^{i-1}) \times p^{k-i}

             =pk+k×(pk−pk−1)~~~~~~~~~~~~~=p^k+k \times (p^k-p^{k-1})

             =(k+1)pk−kpk−1~~~~~~~~~~~~~=(k+1)p^k-kp^{k-1}

同理得：f(pk+1)=(k+2)pk+1−(k+1)pkf(p^{k+1})=(k+2)p^{k+1}-(k+1)p^{k}

f(pk+1)=(k+2)pk+1−kpk−1−f(pk)f(p^{k+1})=(k+2)p^{k+1}-kp^{k-1}-f(p^k)

完整代码：

#include <cstdio>
#include <iostream>
using namespace std;
#define Int register int

const int MAXN = 1000000 , MAXN1 = 1000000;
int x , t , n , m , k , prime[ MAXN + 5 ] , phi[ MAXN + 5 ];
long long f1[ MAXN + 5 ] , f2[ MAXN + 5 ];
bool vis[ MAXN + 5 ];

template<typename _T>
void Read( _T &x ) {
	x = 0; int f = 1;
	char s = getchar( );
	for( ; s < '0' || s > '9' ; s = getchar( ) ) f = s == '-' ? -f : f;
	for( ; s >= '0' && s <= '9' ; s = getchar( ) ) x = x * 10 + s - '0';
	x *= f;
}
template<typename _T>
void Write( _T x ) {
	if( x < 0 ) putchar('-') , x = -x;
	if( x >= 10 ) Write( x / 10 );
	putchar( x % 10 + '0' );
} 

void sieve( ) {
	phi[ 1 ] = 1;
	for( Int i = 2 ; i <= MAXN ; i ++ ) {
		if( !vis[ i ] ) {
			prime[ ++ k ] = i;
			phi[ i ] = i - 1;
		}
		for( Int j = 1 ; j <= k && 1ll * i * prime[ j ] <= MAXN ; j ++ ) {
			vis[ i * prime[ j ] ] = 1;
			if( i % prime[ j ] == 0 ) {
				phi[ i * prime[ j ] ] = phi[ i ] * prime[ j ];
				break;
			}
			else
				phi[ i * prime[ j ] ] = phi[ i ] * ( prime[ j ] - 1 );
		}
	}
	for( Int i = 1 ; i <= MAXN ; i ++ )
		f1[ i ] = f1[ i - 1 ] + phi[ i ];
	
	for( int i = 1 ; i <= MAXN1 ; i ++ )
		for( int j = i ; j <= MAXN1 ; j += i )
			f2[ j ] += i * phi[ j / i ];
	for( int i = 1 ; i <= MAXN1 ; i ++ )
		f2[ i ] += f2[ i - 1 ];
}

long long solve1( int n , int m ) {
    int d = min( n , m ); long long Ans = 0;
	for( Int l = 1 , r ; l <= d ; l = r + 1 ) {
		r = min( n / ( n / l ) , m / ( m / l ) );
		Ans = Ans + ( f1[ r ] - f1[ l - 1 ] ) * ( n / l ) * ( m / l );
	}
	return Ans;
}
long long solve2( int n ) {
	return 2 * f2[ n ] - 1ll * n * ( n + 1 ) / 2;
} 

int main( ) {
	sieve( ); 
	Read( x ) , Read( t );
    while( t -- ) {
        Read( n ) , Read( m );
        Write( n != m ? solve1( n , m ) : solve2( n ) ) , putchar( '\n' );
    }
	return 0;
}
# 莫比乌斯反演
CF17E Palisection 上一篇
P3768 简单的数学题 下一篇
0%